#!/bin/bash

cp protprdict.out proprdict-1.out

sed -i "s/\(^Query:       .+\>\)\(  *.*$\)/\1/g" proprdict-1.out
sed -i "s/^Query:       /Query:       CCDDEF;GHIJKL;/g" proprdict-1.out
sed -i "/^#/d" proprdict-1.out
sed -i '/^Internal pipeline statistics summary:/,/^\/\//d' proprdict-1.out

cp proprdict-1.out proprdicttest

sed -i 's/,/，/g' proprdict-1.out #替换原有的，
cp proprdict-1.out firstdele #复制一份用于制作第一次去除 inclusion threshold 值

#对原始数据进行表格化整理，以用于后期去除 inclusion threshold 值
sed -i '/^Query/,/^Domain/{/^Query/!{/^Domain/!d;};}' proprdict-1.out #删除 Domain 和 Query 开头的行之间的行
sed -i '/^ ---/d' proprdict-1.out #删除表格线
sed -i 's/#/@/g' proprdict-1.out #进一步整理
sed -i 's/!//g' proprdict-1.out #删除表格线
sed -i 's/?//g' proprdict-1.out #删除表格线
sed -i '/^	*$/d' proprdict-1.out 
sed -i '/^ *$/d' proprdict-1.out #删除空白行 


sed -i 's/$/,/g' proprdict-1.out #进一步整理
sed -i ':a;N;s/\n//g;ta' proprdict-1.out #删除换行符
sed -i 's/Query/\nQuery/g' proprdict-1.out #进一步整理
sed -i '/^\s*$/d' proprdict-1.out #删除空白行 
awk '{if(length(max)<length()) max=$0}END{print max > "proprdict2"}' proprdict-1.out #查询最长行并输出到新文件
cat proprdict-1.out >> proprdict2 #因为R语言中读取列数是根据前五行的列数来进行，所以将最长行加入第一行，以防止最长行出现折行

mv proprdict2 proprdict-1.out

cp firstdele ap2_name #复制一份用于制作序列名

#按输出顺序提取序列名
sed -i '/.*/{/CCDDEF;GHIJKL.*\>/!d;}' ap2_name #删除除序列外的所有
sed -i 's/^Quer.*CCDDEF/CCDDEF/g' ap2_name #进一步整理
sed -i 's/\[.*\]$//g' ap2_name #进一步整理
sed -i 's/  *//g' ap2_name #进一步整理,获得名称列表

#R语言处理，以获得带有序列名标识的
R -e 'firstdata <- read.table(file="proprdict-1.out", sep=",",header=FALSE,fill=TRUE)
namesdata <- read.table(file="ap2_name",sep=",",header=FALSE,fill=TRUE)
firstdata <- firstdata[-1,] #去除第一行重复行

numnamecol <- ncol(firstdata) #统计数据列数
names_firstdata <- seq(11,by=1,length.out=numnamecol) #生成列名
names(firstdata) <- names_firstdata #添加列名

firstdata[paste0(names(firstdata),1)] <- namesdata #添加相同列数的蛋白名称列表，并将列名命名为原列名1
firstdata <- firstdata[order(names(firstdata))] #排序表格，使蛋白名称加入每行之后

write.table(firstdata,"first_R_out",row.names=FALSE,col.names=FALSE,na="",sep=",")#输出文件' 

#去除 inclusion threshold 值
#整理R语言输出数据
sed -i 's/"//g' first_R_out
sed -i 's/,CCDDEF/CCDDEF/g' first_R_out
sed -i 's/,/\n,/g' first_R_out
sed -i '/^,$/d' first_R_out
sed -i '/,CCDDEF/d' first_R_out
sed -i '/^,\[ok\]/d' first_R_out
sed -i 's/,//g' first_R_out
sed -i 's/^>> />>/g' first_R_out
sed -i 's/  *$//g' first_R_out

cp first_R_out domain_crip

#循环处理，以删除结构域和基因ID之间的多个单词
wc -l first_R_out > panduan.txt #创建原始判断文件
while [ -s panduan.txt ]; #判断文件是否为空
do
	sed -n 's/\(^>>.*\>\)\(  *.*CCDDEF;\)\(GHIJKL.*$\)/\1_CCDDEF;\3/g'p first_R_out > panduan.txt #提取修改的行
	sed -i 's/\(^>>.*\>\)\(  *.*CCDDEF;\)\(GHIJKL.*$\)/\1_CCDDEF;\3/g' first_R_out #进行行修改
done

#处理 inclusion threshold 值参考文件
sed -i '/^ *$/d' firstdele #删除空白行 
sed -i '/^Domain/,/^Query/{/^Query/!d}' firstdele  #进一步整理
sed -i 's/$/,/g' firstdele
sed -i ':a;N;s/\n//g;ta' firstdele
sed -i 's/Query/\nQuery/g' firstdele
sed -i 's/Scores/\nScores/g' firstdele
sed -i 's/  ------ inclusion/\ninclusion/g' firstdele
sed -i '/^Scores/d' firstdele

sed -i ':a;N;s/\n//g;ta' firstdele
sed -i 's/Query/\nQuery/g' firstdele
sed -i '/^ *$/d' firstdele #删除空白行 

awk '{if(length(max)<length()) max=$0}END{print max > "firstdele_help"}' firstdele #查询最长行并输出到新文件
cat firstdele >> firstdele_help #因为R语言中读取列数是根据前五行的列数来进行，所以将最长行加入第一行，以防止最长行出现折行
mv firstdele_help firstdele

#R语言处理，以获得带有序列名标识的删除目录
R -e 'firstdele <- read.table(file="firstdele", sep=",",header=FALSE,fill=TRUE)
namesdata <- read.table(file="ap2_name",sep=",",header=FALSE,fill=TRUE)
firstdele <- firstdele[-1,] #去除第一行重复行

numnamecol <- ncol(firstdele) #统计数据列数
names_firstdele <- seq(11,by=1,length.out=numnamecol) #生成列名
names(firstdele) <- names_firstdele #添加列名

firstdele[paste0(names(firstdele),1)] <- namesdata #添加相同列数的蛋白名称列表，并将列名命名为原列名1
firstdele <- firstdele[order(names(firstdele))] #排序表格，使蛋白名称加入每行之后

write.table(firstdele,"second_R_out",row.names=FALSE,col.names=FALSE,na="",sep=",")#输出文件' 

sed -i 's/"//g' second_R_out
sed -i 's/,CCDDEF/CCDDEF/g' second_R_out
sed -i 's/,/\n/g' second_R_out
sed -i '/^CCDDEF/d' second_R_out
sed -i '/^$/d' second_R_out #删除空白行 

sed -i '/^inclusion/d' second_R_out 
sed -i '/^Query/d' second_R_out 
sed -i 's/   */,/g' second_R_out
cut -d ',' -f 10 second_R_out > chaifen1.txt
cut -d ',' -f 11 second_R_out > chaifen2.txt
sed -i 's/^.*CCDDEF/CCDDEF/g' chaifen2.txt
paste -d "_" chaifen1.txt chaifen2.txt > firdelete

sed -i 's/$/,/g' first_R_out
sed -i '/^Domain/d' first_R_out 
sed -i ':a;N;s/\n//g;ta' first_R_out #删除换行符
sed -i 's/Query/\nQuery/g' first_R_out
sed -i 's/>>/\n>>/g' first_R_out
sed -i '/^$/d' first_R_out #删除空白行 



cp first_R_out first_R_out-2
sed -i '/^$/d' first_R_out-2 #删除空白行 

cat firdelete | while read line;do
	 sed -i "/^>>$line/d" first_R_out
	 done


sed -i 's/,/\n/g' first_R_out
sed -i '/^$/d' first_R_out #删除空白行 
sed -i '/   @/d' first_R_out
sed -i 's/\[\]//g' first_R_out
sed -i 's/\.\]//g' first_R_out
sed -i 's/\[\.//g' first_R_out
sed -i 's/\.\.//g' first_R_out


sed -i 's/Query:  */Query:/g' first_R_out
sed -i 's/  *\[/\[/g' first_R_out
sed -i 's/  */,/g' first_R_out
sed -i 's/Query/,,,,,,,,,,,,Query/g' first_R_out
sed -i 's/>>/,,,,,,,,,,,,>>/g' first_R_out

cut -d ',' -f 5,9,10,13 first_R_out > secdata
sed -i '1 i ABCD,' secdata #在第一行加入临时辅助标识符
sed -i 's/$/@/g' secdata
sed -i ':a;N;s/\n//g;ta' secdata
sed -i 's/,,>>/\n,,>>/g' secdata

cp secdata secdanames
sed -i 's/_CCDDEF.*$//g' secdanames
sed -i 's/,@,,.*$//g' secdanames
sed -i 's/,,>>/!!/g' secdanames

awk '{if(length(max)<length()) max=$0}END{print max > "secdata_help"}' secdata #查询最长行并输出到新文件
cat secdata >> secdata_help #因为R语言中读取列数是根据前五行的列数来进行，所以将最长行加入第一行，以防止最长行出现折行
mv secdata_help secdata

#提取结构域名称
#R语言处理，以获得带有结构域名标识的删除目录
R -e 'secdata <- read.table(file="secdata", sep="@",header=FALSE,fill=TRUE)
secdanames <- read.table(file="secdanames",sep=",",header=FALSE,fill=TRUE)
secdata <- secdata[-1,] #去除第一行重复行

numnamecol <- ncol(secdata) #统计数据列数
names_secdata <- seq(11,by=1,length.out=numnamecol) #生成列名
names(secdata) <- names_secdata #添加列名

secdata[paste0(names(secdata),1)] <- secdanames #添加相同列数的蛋白名称列表，并将列名命名为原列名1
secdata <- secdata[order(names(secdata))] #排序表格，使蛋白名称加入每行之后

write.table(secdata,"third_R_out",row.names=FALSE,col.names=FALSE,na="",sep="@")#输出文件'


sed -i 's/"//g' third_R_out
sed -i 's/@!!/!!/g' third_R_out
sed -i 's/@/\n/g' third_R_out
sed -i 's/Query/\nQuery/g' third_R_out
sed -i '/^ABCD/d' third_R_out
sed -i '/^,\?!!.*\>$/d' third_R_out
sed -i '/^$/d' third_R_out #删除空白行 
sed -i '/^,*$/d' third_R_out
sed -i 's/!!/,/g' third_R_out
sed -i '/^,,>>/d' third_R_out

sed -i 's/^/DOMAIN,/g' third_R_out
sed -i 's/^DOMAIN,Query:/CHAIN,/g' third_R_out
sed -i 's/\[L=/,1,/g' third_R_out
sed -i 's/\]/,/g' third_R_out

#整理主链
sed -i 's/\(^CHAIN,\)\(CCDDEF.*,\)\(1,.*CCDDEF;\)\(GHIJKL.*$\)/\10.000666666,\3@\2\2/g' third_R_out
sed -i 's/,CCDDEF;@/,/g' third_R_out
sed -i 's/,$//g' third_R_out

#调整序列顺序，并在最前面一列添加标识符，用于后期列名成的生成
cut -d ',' -f 5 third_R_out > firtidyup-2
sed -i 's/^.*CCDDEF/CCDDEF/g' firtidyup-2
cut -d ',' -f 1,2,3,4,6 third_R_out > firtidyup-1
paste -d "," firtidyup-1 firtidyup-2 > firtidyup

sed -i '1 i type,c_EV,begin,end,entryName,order' firtidyup #在第一行加入列名

#制做结构域描述对应文件
cut -d ',' -f 5 domain_crip > domdescrip
sed -i '/.*/{/^>>/!d}' domdescrip
sed -i 's/>>//g' domdescrip


sed -i 's/  /,/g' domdescrip
sed -i 's/CCDDEF.*$//g' domdescrip

sed -i 's/CCDDEF;GHIJKL;//g' domdescrip ap2_name firtidyup

#用R语言进行进一步清洗
R -e 'if (!requireNamespace("BiocManager", quietly = TRUE))
install.packages("BiocManager")

if (!requireNamespace("drawProteins", quietly = TRUE))
BiocManager::install("drawProteins")

if (!requireNamespace("dplyr", quietly = TRUE))
install.packages("dplyr")

if (!requireNamespace("ggplot2", quietly = TRUE))
install.packages("ggplot2")

library(drawProteins)
library(dplyr)
library(ggplot2)

yuandata <- read.table(file="firtidyup",sep=",",header=TRUE)
domdescrip <- read.table(file="domdescrip",sep=",",header=FALSE)
namesdata <- read.table(file="ap2_name",sep=",",header=FALSE)

yuandata <- subset(yuandata,c_EV<=0.01)
yuandata <- yuandata[,-2]

yuandata <- yuandata[,c(1,4,2,3,3,5,4,5,5)] #数据重排
yuandata[,5] <- yuandata[,4] - yuandata[,3]+1 #加入序列长度信息

namesdata2 <- namesdata[,c(1,1)]  #制作含两列序列号的数据



names(namesdata2) <- c("a","b") #添加统一列名
names(domdescrip) <- c("a","b") #添加统一列名

dscriptionall <- rbind(namesdata2,domdescrip) #将序列名称及描述和结构域名称及描述按列（上下连接）合并
dscriptionall <- dscriptionall[!duplicated(dscriptionall),]    #去重

dmolist <- yuandata[,c(2,2)] #提取序列和结构域排序列表


names(dmolist) <- c("dom","crip") #添加统一列名
names(dscriptionall) <- c("dom","crip") #添加统一列名
 #添加统一列名
dscription <- left_join(dmolist,dscriptionall,by = "dom")  #将描述信息对应序列和结构域顺序

names(yuandata) <- c("type","description","begin","end","length","accession","entryName","taxid","order") #添加可用于制图的列名
names(dscription) <- c("entryName","entryName","description") #添加统一列名

decription_zhongzhuan <- dscription[,3] #将描述信息加入数据文件第三列
yuandata[,2] <- decription_zhongzhuan #将描述信息加入数据文件第三列

numnamecol <- nrow(namesdata) #统计数据列数，用于下方制作order列（最后一列，为序列编号）
chain_num <- data.frame(namesdata,seq(1,by=1,length.out=numnamecol)) #生成列名

chain_amio_num <- yuandata[,c(9,9)] #提取order列，用于转换为数字型编号

names(chain_amio_num) <- c("dom","nums") #添加统一列名
names(chain_num) <- c("dom","order") #添加统一列名
data_order <- left_join(chain_amio_num,chain_num,by="dom") #将描述信息对应序列和结构域顺序

order_final <- data_order[3]

yuandata[,9] <- order_final 

final_data <- yuandata

#绘图
p <- draw_canvas(final_data)
p <- draw_chains(p, final_data)
p <- draw_domains(p, final_data)
p <- draw_regions(p, final_data)
p <- draw_motif(p, final_data)
p <- draw_phospho(p, final_data, size = 8) 



p <- p + theme_bw(base_size = 10) + theme(panel.grid.minor=element_blank(), panel.grid.major=element_blank()) + theme(axis.ticks = element_blank(),  axis.text.y = element_blank()) + theme(panel.border = element_blank()) + geom_text(size=800)
pdf("prodomain.pdf",166,6)
p
dev.off()
'

#第一阶段画总图结束，删除过程文件
rm ap2_name chaifen* dom* fir* panduan.txt proprdict* sec* third_R_out